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ABSTRACT 

Low Mach number, high beta fast mode shocks can occur in the magnetic reconnection outflows of 
solar flares. These shocks, which occur above flare loop tops, may provide the electron energization 
responsible for some of the observed hard X-rays and contemporaneous radio emission. Here we present 
new 2D particle- in-cell simulations of low Mach number/high beta quasi-perpendicular shocks. The 
simulations show that electrons above a certain energy threshold experience shock-drift-acceleration. 
The transition energy between the thermal and non-thermal spectrum and the spectral index from 
the simulations are consistent with some of the X-ray spectra from RHESSI in the energy regime, 
E < 40 ~ 100 kcV. Plasma instabilities associated with the shock structure such as the modified-two- 
stream and the electron whistler/mirror instabilities are examined and compared with the numerical 
solutions of the kinetic dispersion relations. 

Subject headings: PIC simulation, solar flares, low Mach number shocks, shock drift acceleration, 
modified two-stream, whistler and mirror instabilities 



1. INTRODUCTION 

Low Mach number (M), high plasma beta (/3 P ) fast 
mode shocks can occur in the magnetic reconnection 
outflows of solar flares. Hard X-ray data from Yohkoh 
and RHESSI has revealed that electrons are energized 
abov e flare loop tops and foot points (e.g., ILin et ahl 
2003). Solar flares are diverse and the associated re- 
connection events are likely "acceleration environments" 
with potentially different mechanisms of particle acceler- 
ation operating on different scales. One potential source 
of particle acceleration, seen from analytic predictions 
(|Blackman fe Field! fl99l and numerical simulations of 
reconnection configurations in which an obstacle was 
present (e.g., iForbesI 119981 iWorkman et al.l 1201 If ), in- 
volves the presence of low Mach number fast shocks in 
reconnection outflows. These "termination" shocks may 
contribute to the high energy acceleration observed over 
some frequency range, studying the potential ways in 
which these shocks can accelerate particles is well moti- 
vated; 

iMann et al.l (p006L l2009h and iWarmuth et"all ({20091 ) 
suggested an electron energization via shock drift ac- 
celerat ion (hereafter SPA) in a sola r flare termination 
shock. IGuo fe Giacalonel (|2010l [20L1 performed hybrid 
simulations for a solar flare termination shock and test 
electrons are effectively energized via the interaction with 
pre-e xisting large-scale magnetic fluctuations. In recent 
work (Park et al. 2012), we performed a full particle-in- 
cell (PIC) simulation for purely perpendicular (i.e. mag- 
netic field perpendicular to shock normal), low M/high 
/3 P shocks, where the upstream magnetic field was di- 
rected out of the simulation plane. We found that both 
electrons and ions are participated in SDA. 

Termination shocks may however deviate from purely 
perpendicularity in solar flares, the extent to which is an 
open question. In the meantime, it is instructive to relax 



the constraint that the shocks are exactly perpendicu- 
lar and consider the case of quasi-perpendicularity. The 
difference is significant because particles can cross back 
upstream for quasi-perpendicular shocks. 

In this paper, we present the results of full PIC simula- 
tions for quasi-perpendicular low M/high j3 p shocks. The 
upstream magnetic field makes angles of 9b = 80°, 82°, 
and 83.5° respectively to the shock normal and the 
shocks satisfy the subluminal condition, V s h/cos8B < 
c, where V s h is the shock speed in the upstream rest 
frame and c is the speed of light. In such a sublumi- 
nal shock, electrons can then be reflected at the shock 
front due to the magnetic mi r ror effect and ga i n en- 
ergy (e.g..lBall fe Meho^l200ll: [Mann et al.l [20061 [2001 
IWarmuth et al.ll2009f K ~ 

One difficulty in performing PIC simulations for SDA 
in solar flares is that they require sufficiently high en- 
ergy electrons above the threshold energy for SDA in 
the simulation to obtain statistically reliable results. 
The commonly used Maxwellian distribution has too 
few above the threshold electrons to be used in a 
simulation. Herein we use a kappa distribution for 
both ions and electrons with k = 10 for the injected 
particles. This increases the number of high energy 
electrons but does not significantly alter the original 
shock structure, which is determined by the bulk of 
the thermal particles. Although the kappa distribu- 
tion is used for computational convenience in this pa- 
per, the kappa distribution can also be physically a 
relevant distribution in the outflow driven electron-ion 
iet (e.g.. lYoon et all 120061: iKasoarova fe Karlickvl 120091: 
IMann et al.ll2006l l2009t IWarmuth et al.ll2009D . 

The goals of this paper are two-fold: (1) studying the 
formation and structure of such shocks, including the tur- 
bulent dissipation mechanism for collisionless shock sus- 
tenance and entropy creation, and (2) studying the par- 
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ticle acceleration mechanism relevant for the soft/hard 
X-ray flux observations in solar flares. We observe the 
same modified-two-stream instability in the shock tran- 
sition region as in the perpendicular shocks (Park et al. 
2012) which can provide the turbulent dissipation in the 
downstream. Furthermore, a temperature anisotropy af- 
ter the shock transition region is found to drive the elec- 
tron whistler/mirror instabilities. 

As we will see, more abundant SDA-accelerated elec- 
trons are found in the present simulations of the quasi- 
perpendicular case compared to the previously sim ula- 
tions of the perpendicular case (jPark et all I2012T ). In 
this context , we a l so ext e nd the theoretical ana lysis in 
IMann et al.1 d20M f2009h : IWarmuth etafl (12001) ' to in- 
clude the electric pot ential jump at the shock front (e.g., 
iBall fc Me lrose 2001) and generalize the electron energy 
spectrum. Both our theoretical analysis and our simu- 
lations show a transition energy, -Etrans.pi between the 
thermal and non-thermal photon spectrum that is de- 
termined by the minimum angle of 9b- We find that 
the transition energy -Etrans.p and the spectral index 8 
from theory and our simulation are consistent with some 
of the X-ray spectra of solar flares from RHESSI (e.g., 
lAltvntsev etaU 120121 ) in the energy range, E < 40 ~ 
lOOkeV. Beyond E = 40 ~ lOOkeV, to maintain the 
power-law distribution up to E ~MeV, additional mech- 
anisms beyond the shock acceleration studied are re- 
quired. 

The rest of the paper is organized as follows. The sim- 
ulation setup is described in section [5J The shock struc- 
ture and particle acceleration are described in section [3J 
We summarize in section 3) 

2. SIMULATION SETUP 

We use the full y relativistic PIC code OSIRIS 
(|Fonseca et al.ll2002l ) to study the formation of and par- 
ticle energization in low-M, high-/3 p shocks, where M 
is the Mach number and j3 p is the ratio of thermal 
to magnetic pressure. To launch a shock, we use the 
moving wall bou ndary condition (jLangdon et al.l 119881 : 
iPark et aLll2012[ ) at the right boundary of the 2D simu- 
lation box. The moving wall method generates a slowly 
propagating shock compared to the more standard fixed 
reflection boundary method, and allows for smaller box 
sizes and more efficient use of simulation time. 

We adopt parameters typical of those fo u nd in 
solar flare reconne ction outflows (jTsunetal 119961 : 
iWorkman et al.l 1201 lh as the upstream conditions for 
our shock. In particular, we use a plasma density 
n = 5 x 10 9 /cm 3 , electron and ion temperatures T e = 
Tf = 0.8kcV(= 9.27 x 10 6 K), and the upstream magnetic 
field strength B = 6G with /3 p = &Tm(T e +T z ) / B 2 = 8.93. 
The magnetic field is in the x-y plane and has an angle 
of 6 B = 80°, 82°, and 83.5° from the shock normals- 
axis) in each simulation. A reduced ion/electron mass 
ratio of rrii/m e = 30 is used to reduce computational 
demands. The Alfven Mach number is chosen to be 
Ma = V\\J Amrriinj B = 6.62, which equates to an up- 
stream plasma flow velocity in the shock rest frame, 
V\ = 0.032c, for mj/m e = 30. If the real mass ratio, 
mi/m e = 1836, is used for the same Ma, then the up- 
stream flow velocity would be 1226km/s(= 0.0041c). The 
super-fast-magnetosonic Mach number M satisfies M = 
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Fig. 1. — The initial energy distribution in the upstream rest 
frame from the simulation. The dashed line is a kappa distribution 
with re = 10 and T = 0.8keV. The dotted line is a Maxwellian 
distribution with T = 0.8keV. 



M A /\J1 + (5/6)/3 p = 2.28. The ratio of the electron 
cyclotron frequency to the electron plasma frequency is 
tt ce /u! pe = 0.0265. With these upstre am values of Ma 
and / 3 P , the Rankine-Hugoniot relation (jTidman fc Kralll 
|1971[ ) for the shocks with the angles used gives the com- 
pression ratio between (2.06,2.50), where the lower and 
the upper limits are calculated for 2D and 3D, respec- 
tively. 

The simulation box is initialized with a kappa- 
distributed ion-electron plasma drifting with Vd = 
0.0213c and T e = T, = 0.8keV, where V d is set to a 
smaller value than the upstream speed V\ (= 0.032c) in 
the shock rest frame considering the shock to be travel- 
ing to the left. In the upstream rest frame, the energy 
distribution is a kappa distribution such as 



f(E) = 
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where T is the gamma function and f(E) is normalized to 
1, L dEf(E) = 1. The kappa distribution approaches 
the Maxwellian distribution as k goes to oo. Figure [1] 
shows the initial energy distribution in the upstream rest 
frame from the simulation. The dashed line is a kappa 
distribution in Equation (JXJ) with k = 10 and T — 0.8keV 
and the dotted line is a Maxwellian distribution with the 
same temperature. The bulk parts (E < 5 keV) of the 
two distributions are very similar. The implementation 
of the kappa distribution in OSIRIS is described in Ap- 
pendix [X] 

A uniform external E z field is set up along the z-axis 
with E x — —VdBy/c. A new plasma of the same kappa 
distribution is constantly injected from the left boundary 
(x = 0) throughout the simulation. The simulation box 
sizes are L x x L y = 250c/o;p e x 50c/u;p e . The grid size 
used is dx = dy = 0.08c/w pe and the time step used is 
dt = 0.056/wp e . For each particle species, 196 particles 
per cell are used. The particle number/cell is fairly large 
to reduce numerical collisions and maintain the kappa 
distribution for a sufficiently long time. 

A linear current deposition scheme is used for all sim- 
ulations in this paper. A periodic boundary condition 
is used in the y-direction for both particles and fields. 
For fields, an open boundary condition is used in the x- 
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Fig. 2. — The ratios, n/ni , Vi/Vx, and B y /B y i , momentum distribution of p x — x, p y — x, and p z —x, y-averaged flow velocity of Vx, V y , 
and Vz, and temperature Tu, T±, T z for the electrons (left column) and the ions(right column) at t = 11200/a>p e . (a) ~ (6) are calculated 
in the shock-rest frame and the other plots are obtained in the simulation frame, (m-o) show the B x , B y , and B z fields in the x — y space. 



direction. Particles that reach x = are re-injected into 
the box with the initial kappa distribution. A t x = L Xl a 
moving wall boundary condition is adopted ( Par k et al.l 
I2012T) and the moving wall speed is set to V^ii = 0.004c. 

3. RESULTS AND ANALYSIS 
3.1. Shock structure 

In Figure [^a&b), we plot the ratios of downstream 
to upstream density n(x)/n\, upstream to downstream 
flow speeds, V\/V{x), and downstream to upstream y- 
magnetic field, B y (x)/B yl . In Figure Hfc-1), we plot the 
phase-space distributions of p x — x, p y — x, and p z — x, 
the y-averaged flow velocity profiles of V x (x), V y {x), and 
V z (x), and the temperature profiles of Tji (x), T±(x), and 
T z (x) (where _L and || are perpendicular and parallel to 
B(x), respectively), for the electrons (the left column) 
and the ions (the right column) at t = 11200/cj pe for the 
9 B = 80° case. 

In this subsection, we present results for the shock 
structure for 9 b = 80° only but other simulations for 
9b = 82° and 83.5° also show similar shock structures. 
The electron energization via SDA shows different energy 
spectra for the different angles of 9b as we will see in the 
next subsection. 

In Figure [2] the shock front is at x ~ lA0c/uj pe and 



moves to the left with a speed of 0.01c in the simula- 
tion frame. In Figure [2ja&b) , the compression ratio is 
r = 2.5 near the shock front and oscillates around 2.0 
downstream, in re asonable agreement wi th the Rankinc- 
Hugoniot relation (jTidman fc KralMl971f) . In Figure [5Jc- 
h) , we show that the electrons and ions are heated down- 
stream, but the ions are conspicuously less heated in the 
y-direction [Figure [2f]. In Figure HJi&j), the flow ve- 
locity, V x = 0.021c, in the upstream decreases down- 
stream to the moving wall speed, V wa \\ = 0.004c. The z- 
component of the velocity, V z , is oscillating around zero 
downstream due to the x B y -drift while V y ~ 0. In 
Figure [5Jk&l), we see that the downstream temperature 
is anisotropic, that is Tj_ > Tji with T e j_ — T ez 1.66keV 
and T e ii s» 1.5keV for the electrons, and Ti± = Ti Z s» 
2.5keV and Tq w l.OkeV for the ions. 

Figure [^m-o) show B in x — y space. The B x and 
B z fields show oscillatory patterns along the ?/-axis near 
the shock front (x ~ 150c/o;p e ) and the B y field shows 
a rippled surface along the y-axis. These variations are 
due to electron whistler/mirror instabilities driven by the 
electron temperature anisotropy, T e ± > T e u. 

Figure [3] shows the ^-averaged E, B, and the potential 
energy |e|$(x)(= e J E x dx) at t = 11200/w pe , mea- 
sured in the shock rest frame. The E y , E z , and B x fields 
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Fig. 3. — jy-averagcd E and B fields, and the potential energy \e\<f>(x) at t = 11200/wpe measured in the shock rest frame. 
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Fig. 4. — (a)~(d): Fourier spectra of the E and B fields in the k x — k y space from the simulation. The Bo field has an angle of 86° from 
the x-axis. (e&f): The growth rates for the electron whistler /mirror modes from the kinetic dispersion relation (Appendix [Cj , (f): The 
growth rates for the electron whistler/ion cyclotron modes for k X Bq = 0. 



are approximately constant across the shock while E x 
and B z oscillate with a wave number k = 0.2cj pe /c down- 
stream. The potential energy jump |e|A$ at the shock 
front is ~ 3.5keV [Figure [3Ji] . Some ions in the low en- 
ergy tail reflect at the shock front and drive the modified 
two-stream instability. 

In Figure SJa-d), we plot the Fourier spectra of the E±, 
B±, Bu, and B z fields, where _L and || are perpendicular 
and parallel to Bo, respectively. The Bn(= Bq x x + Bo y y) 
is an averaged field over 150 < x < 170(c/u! pe ) in 
Figure [5] and has an angle 86° from the £-axis. The 
electron whistler and mirror modes drive n by electron 
temperature anisotropy, T e ± > T e u (e.g.. iGarvl 119931 : 
IGarv fc Karimaba di 2006) are indicated. The dominant 
electron whistler mode is along the Bo with k = 0.5ui pe /c. 
The dominant mirror mode has k — 0.8u pe /c and makes 



an angle 6 = 21° from Bo- The electron whistler modes 
are circularly polarized and therefore not seen in the Bu 
plot [Figure Ht] . 

To compare the simulation results with the kinetic dis- 
persion relations, we numerically solve the dispersion re- 
lation (Appendix |C|) in Figure HJe&f). Here we assume 
that the electron and ion distributions are bi-Maxwellian 
for simplicity. We use the parameters extracted in the 
region of 150 < x < 17Q{c/uj pe ) at t = 11200/^ in Fig- 
ure CJk&l) with T e u = 1.33keV and T el _ = 1.8keV for 
the electrons, and Tjii = 0.85keV and Ti± — 2.3keV for 
the ions. The magnetic field strength Bq — 16G. The 
maximum growth rate for the electron whistler mode 
is u>i = 5.2 x 10 _4 w pe at fcii = 0.47(w pe /c) [Figure 
d^&f ] . The maximum growth rate for the mirror mode is 
LOi = 1.4 x 10~ 3 w pe at k± = 0.63 and ku = 0.8&(uj pe /c). 
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Fig. 5. — (a&b): Ion and electron distributions in the shock 
transition region, 155 < x < 175(c/u;p e ) at t = 8400/<^p e - We 
fit the distributions with Maxwellian distributions (dashed lines), 
(c&d): The growth rates and real frequencies obtained from the 
MTSI dispersion relation (Appendix [BJ . 



The solutions for the dispersion relation when k x B = 
(Appendix [D]) and the ion-cyclotron instability weaker 
than the electron whistler instability are explicitly shown 
in Figure HJf). 

Figure [Sf a&b) shows the electron and ion distributions 
in the shock transition region in 155 < x < 175(c/w pe ) 
at t = 8400/cjp e . We observe that 27% of the incom- 
ing ions are reflected at the shock front [Figure [5^] and 
the modified two-stream instability (MTSI) can be ex- 
cited. The MTSI is electrostatic with its wave vector 
along the cc-axis(Appendix[B]). The electrons have a tem- 
perature of T ex — 1.65keV and drift with V ex = 0.0064c 
in the simulation frame [Figure [5})]. We fit the dis- 
tributions using drifting Maxwellians. The electrons 
are magnetized with B — 10G while the ions are as- 
sumed to be non-magnetized. In the electron rest frame, 
the drift velocities for the incoming and reflecting ions 
are V x \ n = 0.0116c and V xie = —0.0134c, respectively. 
Both incoming and reflecting ions have a temperature of 
T xin = T XTe = 0.98keV. In Figure Efc&d) , we solve the 
kinetic dispersion relation for the MTSI. The maximum 
growth rate is uii = 3.2 x 10~ 4 w pe at k = 0.2c/cj pe and 
the real frequency is u r — — 0.001 5w pe . 

Turbulent dissipation is needed to random ize the up - 
stream flow to form collisionless shocks fe.g.. lWuH l982). 
The macroscopic jump conditions across a shock are 
essentially independent of the source of microphysical 
turbulence, as long as there is such a source. Plasma 
instabilities in the shock transition region are natural 
sources for this needed turbulent dissipation and we 
consider three possible instabilities: (l)lower hybrid in- 
stability (the excited modes a re in the B x VB di- 
rection) (e.g.. iZhou et al.l I1983T). (2)whistler/mirror in- 
stability, and (3)M TSI (e.g.. iPapadopoulos et all 119711 : 
iWagner et al.lll97l"l ). The lower hybrid instability is pre- 
cluded in these 2D simulations where the B-field is in 
the x — y simulation plane. In addition, we notice that 
the shock structure (i.e. the compression ratio) in the 



simulation with an in-plane B field is the same as that 
with the out-of-plane B field which then precludes the 
whistler/mirror instability. We are led to conclude that 
the MTSI is the most likely candidate to provide the 
needed turbulence fe.g. JPark et al.ll2012( l. 

3.2. particle heating via shock drift acceleration 

In Figure [6j we plot the electron and ion energy distri- 
butions, f(E)(= dN/dE) upstream(50 < x < 155c/w pe ) 
and downstream(155 < x < 250c/w pe ) at t = 8400/w pe 
for 9 B = 80°. We fit the thermal (bulk) part of the 
distributions using a kappa distribution with k = 10 
(dashed lines), yielding T e = Tj = 0.8keV upstream and 
T e = 1.58keV and Tj = 1.8keV downstream. In Figure 
[UJa&d), we plot the electron and ion phase-space dis- 
tributions in E — x, overlaid with B y field to indicate 
the location of the shock front. Abundant non-thermal 
electrons are seen ahead of the shock at x = 155c/w pe , 
traveling as far back asi= 10c/uj pe [Figure HJi]. 

In Figure [6th), the electron energy spectrum shows 
a deviation from a thermal distribution at E ~ 12kcV 
and the spectrum becomes steeper at E ~ 50keV. The 
dotted line in Figure [6^b) is the theoretical energy dis- 
tribution via SDA when the electric potential energy 
jump |e|A$ w 3.5(kcV) at the shock front is consid- 
ered (See section I5.2.1j) . In Figure [SJe) , the ion energy 
spectrum shows a deviation from a thermal distribution 
at E ~ 15keV and the spectrum becomes steeper at 
E - 30keV. 

In Figure [TJa-d) , we plot a typical track of an electron 
experiencing SDA, overlaid with the shock front and the 
subsequent compression peaks [Figure [7^]. The electron 
is reflected at the shock front [Figure^] and gains energy 
from 8keV to 50keV [Figure Eh]- After the reflection, the 
electron drifts along the upstream magnetic field lines 
[Figure [7b] ■ In Figure E{d), we plot the electron in the 
vj_ — vu phase-space. The electrons in the region I are 
transferred to the region II after the reflection. 

In Figure [TJe&f), we plot a typical track for an ion 
experiencing SDA. Whether or not an ion gains energy 
via SDA depends on its inc ident speed and angle of in- 
cidence at the shock front ([Kirk I Il994h . When the ion 
meets the shock front, it turns back toward the upstream 
with a larger gyro-radius [Figure [7^,] and is accelerated 
by the E z field. The kinetic energy of the ion increases 
from 5keV up to 20keV [Figure [7b] ■ A detailed analysis 
for ions experiencing SDA in p erpendicular sho cks was 
described in our previous work ([Park et alJl2012ft . 

3.2.1. Electron spectrum via SDA 

In this subsection, we g eneralize the electron en - 
ergy spectrum via SD A by iMann et al.l (|2006l 12009( 1: 
iWarmuth et al.l ([2009( 1 to include the el ectric potential 
jump , eA$, at the shock front (e.g., I Ball fc Melrose! 
2001) and compare with the PIC simulation results. 

First, we consider the de Hoffmann- Teller(dHT) frame 
(denoted by ') where the motional E z (= —Vi/cB y ) field 
vanishes. The dHT frame is obtained by boosting with 
v s = V s h/cos0B along the magnetic field line in the up- 
stream rest frame. (Here we consider the negative shock 
speed, V s h < 0, for a shock traveling to the -x direction.) 
The maximum 8b for the existence of the dHT frame is 
given by 6> Bmax = cos~ 1 (l/ sh /c)(< 90°). 
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Fig. 6. — (a&d): Energy distribution vs. rc-ranges with the B y field is over-plotted with an arbitrary scale. The electron(b&c) and 
ion(e&f) energy distributions upstream (50 < x < 155c/uj pe ) and downstream (155 < x < 250c/u pe ) at t = 8400/cjpe- In (b), the dotted 
line is a theoretical energy distribution when the potential energy at the shock front is e<l> = — 3.5(keV). We fit the thermal distributions 
with kappa distributions with k = 10(dashed lines). 
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Fig. 8. — The condition for SDA in the 0» — /3j_ space 
when p s = —0.293. The incoming electrons in the region 
I(light sh aded) are transfer r ed to the region II(dark shaded) after 
reflection! Mann ct al. 2006, 2009). The hyperbolic curve is given 
by the potential energy e3> = — 3.5keV at the shock front, which 
approaches the straight dotted lines when e<J> — ► 0. 



In the dHT frame, the conditi on for an electron to 
reflect at the shock fro nt is (e.g., iBall fc Melrosd I200H: 
IMann et alJl2"006l [2009) 



f3' ± > /3||tana and P'n > 0, 



(2) 



where _L and || are perpendicular and parallel to the 
upstream Bi, respectively, j3',, ± = v',, j_/c, and ao = 



sin - -JB1/B2 (~ siri^-i/l/r). If we include the po- 
tential energy e<I>(< 0) at the shock front, the reflec- 
tion condition, Equation ( [2t, can be written in the non- 
relativistic limit by (e.g.. IBall fc MelrosdfeOOll ) 




2e$ 
m e c 2 J 



jtanao, /3[| > 0. 



(3) 



For relativistic electrons with /3m ^> ^2\e\<& /m e c 2 (= 

0.117), the potential energy e$ is negligible and Equation 
((3|) is reduced to Equation ((2]). 

Using the Lorentz transformation between the dHT 
frame(') and the upstream rest frame, 

011 ~ & a, Pi 



P ' 1 ^ 7.(1 -V-)' 



(4) 



where f3 s = v s /c(= V sh / {ccos8 B )) and 7 S = - 
Equation is transformed in the upstream rest frame 
into 



/3_l > 7 s tanao 
P\\>Ps- 



(P\\ - Ps 



2e$ 



1/2 



(5) 



The velocity of the electron after reflection in the dHT 
frame is given by /3'., = — (3L and /3' rX — fi' i± , where 
the indices i and r represent the incoming and reflecting 
elect ron, respectively. In the upstream rest frame, one 
gets (jMann et allHOOl 20Q9|) 



Pi 



2ft-/3 r ||(l + /3 2 ) 
1-2^0.+ ft 



Pi. 



PrJ 



7s 2 (l-2/3 r ||/3 s +/?2 



(6) 
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Fig. 9. — The transition energy -Etrans.e vs 8b for e<t> = —3.5 
keV(solid) and keV (dashed) from Eq.(|8]l. The energy point i?2,e 
vs G B (dot-dashed) from Eq.©. Here V sh = 0.0041c, r = 2.5 and 
mi/m e = 1836 (M A = 6.62 and P = 8.93.) 



In Figure [FJ we plot the SDA condition in the j3\\ — j3± 
space for (L = —0.293 as an example (e.g., IMann et al.l 
120061 12009ft . Equation (J5J defines the region I(light 
shaded region) and the hyperbolic curve is given by the 
potential energy e$ = — 3.5keV at the shock front. When 
e<& goes to zero, the hyperbolic curve approac hes the 
straight lines as seen in IMann et al.1 ([2006, 2009]). Equa- 
tion ([6]) implies that the incoming electrons in the region 
I are transferred to the region II (dark shaded region) af- 
ter reflection. As 9b increases to 6>B max (= cos _1 (V r s / l /c)), 
f3 s goes to —1 in Figure[8jand the number of electrons sat- 
isfying the reflection condition in Equation ([5]) decrease 
to zero. Therefore, no electron reflects at the shock front 
in the super luminal shocks (where \(3 S \ > 1) or perpendic- 
ular shocks. 

The threshold energy -Ethrcs f° r SDA is given by the 
shortest distance from the origin to the hyperbolic curve, 
oq(= Ahrcs) in Figure [SI namely 



P- 



thr 




P 2 



thrcs 



1 m e c , where 



(7) 



D 2 [(1 - Ps) 2 D 2 P 2 {P - (1 + P 2 )} + PI - P] 
l + D 2 (l-P/3 2 ) 



where D = 7 s tan«o and P = 2e$/m e c 2 . 

The transition energy £"trans,e between thermal and 
non-thermal electron populations is given by the distance 
from the origin to the point p, op(= /3trans), in Figure [8J 
such as 



E, 



Pt 



1// 
Pi - (1 



1 



B 2 - 1 

A^trans 



m e c 



where 



(8) 



/3 2 )(2e$/m e c 2 )tan 2 a . 



The energy point £?2,e where the maximum energy ra- 
tio of t he reflecting to t he in coming electron occurs is 
given bv lBall fc Melrosel (pOOl in the non-relativistic ap- 
proximation, namely 



E, 



2,e 



PI 



1 m e c 



B 2 = 2/3 s cos(a /2), (9) 

Beyond E = i?2,e, the spectral index <5, defined as in 
f(E) oc E~ s , becomes increasing. Figure [9] shows how 
the transition energy E trans ,e for e$ = —3.5 keV(solid) 
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Fig. 10. — (a-c): The upstream electron distribution in the /3m — 0± plane from the simulations for 9b = 80°, 82°, and 83.5°. (d~f): 
The upstream electron energy distributions in the upstream rest frame from the theoretical distributions with e$ = O(dashed) and 
e<E> = — 3.5keV(dot-dashed), and the simulations (solid). We fit the thermal distribution with a kappa distribution with T = 0.8keV and 
ft = 10(dotted). Here V ah = -0.032c, r = 2.5 and m,i/m e = 30 (M A = 6.62 and /3p = 8.93). 



and keV(dashed), and E% e (dot-dashed) varies with 9 b 
when V sh = 0.0041c, r = 2.5 and m,/m e = 1836 [M A = 
6.62 and /3 P = 8.93). 

The reflected electron distribution in the upstream rest 
frame is written as 



/ r (/3 r ) = / i [/3 i ( / 3 r )](d 3 A/rf 3 ^r)e(/3 s - /3 r |,)x 



© /3r_L - 7 s tana 



2e$ 

m e c 2 



l/2\ 



(10) 



where is the step function and the term d 3 /3i/d 3 (3 r is 
given by the Jacobian determinant, 



d 3 /3 t _ A 

d 3 p r PrX 



d(Ax,ft||) 



(l-2/3 r[i /3 s +/3 s 2 ) 4 - 



(11) 



Here we let the incoming distribution i n Equation (|10l) be 
a sem i-relativistic kappa distribution (jMann et al. 2006, 
2009), 



/,(&) = c R 1 + 



2(7—1) m e c 2 
2k- 3 T~ 



(12) 



where 1% = 1/^1- 0\ and c K = 1/ / d 3 A/i(ft)- 

We calculate the upstream electron energy distribu- 
tions, /(-E) = f((3)d(3/dE, in the upstream rest frame 



using Equations (|6"l fT2")) to obtain 



for /3 < /Sthros 

for /3 t hrcs < j9 < A r, 



for /3 > A runs 



(13) 



where i = cosf5 and 9 is the electron's pitch angle-the 
angle between (3 and Bi. The boundaries ti(= cos#i), 
£a(= cosf5 2 ), and t3(= COS6I3) as shown in Figure [8] are 
given as t% = /3 s /fl, and t\ and t§(< t±) are roots of the 
equation, 

t?, 3 {1 + £ 2 (1 - P/? s 2 ) } - 2^ 2 § (1 - P) 







+^2-(4 2 -^)-l=0. (14) 

In Figure [TUta-cl , we plot the upstream electron dis- 
tributions in the /3m — /3j_ plane from the simulations for 
9b = 80, 82 and 83.5°. We see that electrons in region I 
are transferred to region II. As 9b increases, the number 
of electrons which participate in SDA decrease since the 
energy threshold for SDA in Equation (fTJ) increases. 

In Figure [TUJd^f), we plot the normalized upstream 
electron energy distributions, f(E)(= 1/NdN/dE), in 
the upstream rest frame for 9b = 80,82 and 83.5°. 
We compare the results of Equation (Ti~3"]) for e$ = 
O(dashed) and — 3. 5(keV) (dot-dashed) with the simu- 
lation results(solid). The dotted line is the incoming 



Particle-in-cell simulations of particle energization via shock drift acceleration 



9 




£(keV) £(keV) 

Fig. 11. — (afe b) : The photon distribution via bremsstrahlung radiation for 8g = 80°, 82°, and 83.5° when the electron distributions 
are given by Eq, ]13H and by simulations, (c): The averaged photon distribution of 80° < 8b < $Bmax(= 88.17°) from the theoretical 
distributions with e<E> = and — 3.5keV. (d): The averaged photon distribution of 8b = 80°, 82°, and 83.5°. The dotted line is the 
photon distribution when the electron distribution is a kappa distribution with T = 0.8keV and k = 10. Here = —0.032c, r = 2.5 and 
m,/m e = 30 (M A = 6.62 and /3 P = 8.93). 



kappa distribution with T = 0.8keV and k = 10. Us- 
ing Equation ((SJ, the transition energy points are given 
by £ t rans,e = 11.3keV s 16.5keV, and 24.2keV for the 
angles, 9b = 80°, 82°, and 83.5°, respectively, when 
e$ = — 3.5keV. Using Equation ((9j, the energy points 
i?2,e's, beyond which the spectral index increases, are 
given by i?2, e = 34keV, 56keV, and 93keV for the an- 
gles, B = 80°, 82°, and 83.5°, respectively. Here 
V sh = -0.032c, r = 2.5 and mi/m e = 30 (M A = 6.62 
and P P = 8.93). 

3.2.2. Bremsstrahlung radiation 

Electrons accelerated via SDA move along the mag- 
netic field lines as can be seen from the electron track- 
ing in Figure EJb) and will collide with ions to emit 
bremsstrahlung radiation. In the low-frequency limit of 
uib/^v -C 1, where ui is a photon angular frequency, b 
is the impact parameter, v is the electron speed, and 
7 = 1/ yj\ — (v/c) 2 , the number of photons per unit 
frequency per unit volume per un it time produced by 
the b remsstrahlung process is fe.g.. lRvbicki fc Lightmanl 
fl979h 

dN = 16Z 2 e 6 n e n, f c 7^ hm e v 2 \ 
dwdVdt 3c 3 m 2 h 2 uj J v . v \ haj J ' 

where f(v) is a normalized electron distribution and w m i n 
is determined by the equation ftw = (7 m i n — l)mc 2 . 

Given the electron distributions in Figure nUT d-D. we 
calculate the number of photons per unit time(sec) per 
unit energy(keV) using Equation (1151) for the ion den- 
sity of rii — n e and the volume of the region emitting 



the X-rays, V — 10 27 cm 3 . In Figure [TTJa&b) , we show 
the photon distribution for 0b = 80 (solid), 82 (dashed), 
and 83. 5 (dot-dashed) when the electron distribution is 
given by the theoretical distribution in Equation (|13p 
with e<I> = — 3.5keV [Figure ITTa ] and by the simula- 
tion [Figure ITTb]. The dotted line is the photon dis- 
tribution when the electron distribution is given by a 
kappa distribution with T = 0.8keV and k = 10. Here 
V sh = -0.032c, r = 2.5 and m l /m e = 30 {M A = 6.62 
and P P = 8.93). 

In Figure 111( c) , we plot the averaged photon distribu- 
tion of 80° < B < 6 , Bmax(= 88.17°) from the theoretical 
distribution in Equation (Tf"3"|) with e<i> = (dashed) and 
e$ = — 3.5keV(solid). In Figure ITlTd). we plot the aver- 
aged photon distribution of Ob — 80, 82, and 83.5° from 
the electron distribution given by Equation (|13|) with 
e<t> = O(dashed) and e<I> = — 3.5keV(dot-dashed), and 
by the simulation(solid). In Figure fTTT d). we notice that 
the theoretical result with e$ = — 3.5keV(dot-dashed) is 
in good agreement with the simulation result (solid). 

In Figure fTlTd). the simulation result shows that the 
transition energy -Etrans.p between the thermal and non- 
thermal photon spectrum is -Etrans.p ~ lOkeV and the 
energy point E2 lP , beyond which the photon spectrum 
becomes steeper, is E^, p ~ 40keV. The spectral index 
8 is 8 = 3 (simulation) and 8 = 2. 2 (theory with e<& = 
-3.5keV) in 10keV< E < 40keV and 8 = 7(simulation) 
in E > 40keV. For emission from multiple shocks with 
different #b's, the transition energy £"trans,p and E% jP 
would be dominated by the shock with the minimum 0b- 
Note that the transition energy -Etrans.p for the photon 
distribution is a bit smaller than i?trans,e for the elec- 
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Fig. 12. — The averaged photon distribution of 0Bmin l£ ®B < 
0Bm ax(= 89.77°) from the theoretical distributions in Equation 
l|13|l with e<E> = — 3.5keV for several minimum values of 8b- Here 
V sh = 0.0041c, r = 2.5 and mi/m e = 1836 (M A = 6.62 and 
0p = 8.93). The transition energy increases as 9s m i n increases. 



tron distribution because bremsstrahlung photons are 
produced by elec trons with high er energies than the pho- 
ton energy fe.g.. iHolman I [2003). The energy point E^, p 
is approximately given by Equation 

RHESSI data for several solar flares (Table I in 
lAltvntsev et al.l (|2012| )). shows that the spectral index 5 
is in the range 2.5 < 5 < 3 and the transition energy is in 
the range 12.1 keV < -E t rans, P < 29.2 keV. Therefore, the 
electron energization via SDA well explains some of the 
RHESSI X-ray spectra for the energy regime, E < i?2, P 
and how the transition energy is related with the shock 
geometry, i.e., the minimum 0b- 

The observed RHESSI spectra do not show a steepen- 
ing beyond E — i?2,p, and thus the theory herein can- 
not by itself account for the electron acceleration to pro- 
duce those photons. This indicates additional mecha- 
nisms, such as diffusive shock acceleration, are required 
for further electron energization to maintain the power- 
law spectrum up to E ~ MeV. It is not unreasonable 
to expect that solar flares involve multiple acceleration 
mechanisms operating on a range of scales. 

3.3. Effects of a realistic proton/ electron mass ratio on 

the spectra 

When the actual ion/electron mass ratio of mi/m e = 
1836 is used, the shock speed in the upstream rest frame 
is reduced by a/30/1836 compared to our simulations 
when the Mach number M and the plasma (3 P are fixed. 
The compression ratio r and the electric poten tial jump 
eA$ are unchanged for a fixed M and (3 p (e.g., iHoshinol 
2001). Therefore, the shock structure are not expected 
to change for a realistic mass ratio simulation. 

The electron energy spectrum in Equation ([6]) de- 
pends only on j3 s , eta, an d e<fr. For a fixed M and /3 P , 
a (= sin -1 \/ 'Bij B 2 ) and e4> are unchanged. From the 
definition of /3 S (= V s u / ccosO b) , we only need to change 
the angle from 9b into 0' B to make the electron energy 
spectrum unchanged for the mass ratio mi/m e = 1836 



cos9b 
a/1836/30 



(16) 



For example, 80° < B < 6» B max(= 88.17°) for mjm e = 
30 corresponds to 88.73° < 9 B < 89.77° for mi/m e = 
1836. 

In Figure [121 we plot the averaged photon distribution 
of #Bmin < 0b < #Bmax(= 89.77°) for a real ion/electron 
mass ratio, m,i/m e — 1836, from the electron distribu- 
tion given by Equation (|13|) with e<& = — 3.5keV. Here 
V sh = 0.0041c and r = 2.5 (M A = 6.62 and j3 p = 8.93). 
The different values of 6>Bmin, from Bmin = 88.73° to 
89.43°, give the different transition energy points, from 
-Etrans.p = 10 to 35 keV. The power indices of the 
photon spectrum are nearly the same as 5 ~ 2.2 in 
-E^rans p < E < ^2p, where E<x p runs from 40 to 150 
keV. ' 

4. CONCLUSION 

In summary, we studied quasi-perpendicular, low 
M /high j3 p shocks with full PIC 2-D simulations using a 
reduced ion/electron mass ratio mi/m e = 30. The shock 
compression ratio we found was in agreement with the 
Rankine-Hugoniot relation. Electron whistler/mirror in- 
stabilities driven by downstream temperature anisotropy 
were observed. A modified two-stream instability due to 
the incoming and reflecting ions in the shock transition 
region was also observed. 

Abundant non-thermal electrons accelerated via SDA 
were observed upstream. We compared the electron 
energy distribution from the simulation with the dis- 
tributions derived by extending a theoret ical model 
(|Mann et al.ll20 06. 2009; iWarmuth et al]l2009h to include 
electric potential energy at the shock front, and found 
good agreement. 

In perpendicular shocks, however, SDA can be 
achieved only by particles transmitted into the down- 
stream, and the energy obtained by the transmitted par- 
ticles are smaller tha n that by the reflected ones (e.g., 
iBall &: Me lrose 200l|). Therefore, such abundant non- 
thermal electrons observed in quasi-perpendicular shocks 
in this paper are not seen in perpendicular shocks (e.g., 
IPark et al] 12012^ . 

We calculated the photon flux via bremsstrahlung ra- 
diation from the electron distributions from both the the- 
ory and simulation and we identified a transition energy, 
-£<trans,p marking the change between a thermal and non- 
thermal photon spectrum. -Et ranSiP is determined by the 
minimum 9b from multiple shocks with different #b's in 
#Bmin < 0b < 90°. Different solar flares have different 
^Bmin's in their termination shocks and therefore show 
different transition energy points. 

From the simulations, the averaged photon spectrum 
of B = 80°, 82°, and 83.5° gives a spectral index S ~ 3 
in 10keV< E < 40keV and the spectral index increases 
beyond E = 40kcV. The spectral index S ~ 3 as well 
as the transition energy, -Btrans,p = lOkeV, well explains 
some of the RHESSI X-ray spectra in the energy regime, 
E < 40keV. To account the spectral index of the RHESSI 
X-ray spectrum up to E ^MeV, however, additional 
mechanisms other than SDA are required. 



Particle-in-cell simulations of particle energization via shock drift acceleration 



11 



Note that although our simulations were performed us- 
ing m,i/m e = 30, we analytically scaled the results to the 
realistic ion/electron mass ratio, and found that the pre- 
dicted photon spectra are indeed insensitive to this mass 
ratio. 
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APPENDIX 

GENERATION OF A KAPPA DISTRIBUTION 

To generate a kappa distribution, we use the random 
number distribution by iLeitner et al.1 (|2011| ) . We denote 
{Ni(0, = 1,2, N} as a set of normally distributed 
random numbers with the mean of and the deviation 
of 1, and {Ui(0, l)\i — 1,2, N} as a set of uniformly 
distributed random numbers between [0,1]. Then a se- 
quence of the random number Xi given by 



Xi = b 1 N i (0,l) + b 2 U i (0,l)x 



i-i (Al) 

generates a ID kappa distribution with the mean of 
and the deviation of 1, and the c oefficients, 6i and hi 
determine the index k (Figure 8 in ILeitner et atll201lH . 

To implement the kappa generator in OSIRIS, we de- 
termine the particle's initial momentum p as 

Pi = PthXi + Pd, (A2) 

where Xi is given by Equation (jAll) . p t h and pd are the 
thermal and the drift momentum, respectively. For k = 
10, we choose b\ — 0.5 8 and &2 = 1.15. Then a sequence 
of Pi in Equation (|A2j) generates a 3D kappa distribution 
in Figure [T] 

DISPERSION RELATION FOR THE 
MODIFIED TWO-STREAM INSTABILITY 

The kinetic dispersion relation for electrostatic insta- 
bilities is 

1+ K s (k,u)=0, (Bl) 

s— e,i 

where K s (k, oj) is the susceptibility. For magnetized elec- 
trons with i sotropic Ma xwellian distributions, Ae(k, oj) 
is given by (|Garv II1993D 



K e (k,oj) 



oj,, 



J eth' 



m— — oo 



(B2) 

where v e th is the thermal velocity, I m is the modified 
Bessel function of the first kind, A e = k\v 2 th /Q 2 , VL e = 
eB /(m e c)(< 0), Z is the plasma dispersion function, and 

C e n = (oj - mfi e )/( v / 2fcj| v e th)- For unmagnetized ions 
with drifting Maxwellian distributions, 



«i(k,w) 



pi 



(B3) 



where = (oj — k • Vy) / (s/2kvith) and ' is the derivative 
with respect to £ s . 



Here we consider that B sa B y y, k = k x x V ' s d — 
V s dX. Then the term £e-Z(£™ 1 ) in Equation (|B2I) becomes 
—oj/(oj — mfl e ) as fc|| goes to 0. The dispersion relation 
for the MTSI becomes 



V lth k2 



1 - e 



J2 J ™(Ae 



E 



ps 



Z'(6) = 0.(B4) 



DISPERSION RELATION FOR 
BI-MAXWELLIAN DISTRIBUTIONS 

The kinetic dispersion relation for electromagnetic in- 
stabilities is 



det ( hkj - k 2 Sij + ^jdij + k 2 ^ S sij = 0, (CI) 



where S s ij is the conductivity tensor. Here we consider 
that k = k x x + k z z, Bo = Bqz, and the electro ns and 
ions h ave bi-Maxwellian distributions. Following IGarv I 
(|1993f) . we get 



<4 e- A ° 
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(C3) 

E rnI m (X s )Z'(C s 



m— — oo 



-Wi^- 1 ) E m 2 I m (X s )Z<(C) 



(C4) 
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S a 
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(C7) 
(C8) 



where \ s = k\v% h /fl 2 and £™ = (uj-mVl s )/{^j2k\\v sth ). 
When k x Bq = 0, the dispersion relation, Equation 



CTl). reduces to 



(C9) 



where Sjr is given by 



k 2 c 2 



1 / T. 



(CIO) 
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